home *** CD-ROM | disk | FTP | other *** search
/ Aminet 6 / Aminet 6 - June 1995.iso / Aminet / gfx / 3d / irit50src.lha / irit5 / symb_lib / morphing.c < prev    next >
Encoding:
C/C++ Source or Header  |  1995-01-30  |  28.7 KB  |  772 lines
  1. /******************************************************************************
  2. * Morphing.c - A simple too to morph between two compatible surfaces.          *
  3. *******************************************************************************
  4. * Written by Gershon Elber, Jul. 92.                          *
  5. ******************************************************************************/
  6.  
  7. #include "symb_loc.h"
  8. #include "genmat.h"
  9. #include "geomat3d.h"
  10.  
  11. #define LINE_EPSILON 1e-3   /* Tolerance of a curve to be considered a line. */
  12.  
  13. static void InterpolateVector(CagdRType x1,
  14.                   CagdRType y1,
  15.                   CagdRType x2,
  16.                   CagdRType y2,
  17.                   CagdRType *x,
  18.                   CagdRType *y,
  19.                   CagdRType Blend);
  20. static void TangencyFilterCrvMorph(CagdCrvStruct **CrvSeq1,
  21.                    CagdCrvStruct **CrvSeq2);
  22. static CagdCrvStruct *CrvMorphSeqCornerCut(CagdCrvStruct *Crv,
  23.                        CagdRType MinDist);
  24. static CagdCrvStruct *CrvMorphOneCornerCut(CagdCrvStruct *Crv,
  25.                        CagdRType Blend);
  26. static void DistanceFilterCrvMorph(CagdCrvStruct *MorphSeq,
  27.                    CagdRType MinDist);
  28. static CagdRType *ComputeCrvCentroid(CagdCrvStruct *Crv);
  29. static CagdRType EstimateIntDistCrvCrv(CagdCrvStruct *Crv1,
  30.                        CagdCrvStruct *Crv2);
  31.  
  32. /*****************************************************************************
  33. * DESCRIPTION:                                                               M
  34. * Given two compatible curves (See function CagdmakeCrvsCompatible),         M
  35. * computes a convex blend between them according to Blend which must be      M
  36. * between zero and one.                                 M
  37. *   Returned is the new blended curve.                         M
  38. *                                                                            *
  39. * PARAMETERS:                                                                M
  40. *   Crv1, Crv2:  The two curves to blend.                                    M
  41. *   Blend:       A parameter between zero and one                            M
  42. *                                                                            *
  43. * RETURN VALUE:                                                              M
  44. *   CagdCrvStruct *:   Crv2 * Blend + Crv1 * (1 - Blend).                    M
  45. *                                                                            *
  46. * KEYWORDS:                                                                  M
  47. *   SymbTwoCrvsMorphing, morphing                                            M
  48. *****************************************************************************/
  49. CagdCrvStruct *SymbTwoCrvsMorphing(CagdCrvStruct *Crv1,
  50.                    CagdCrvStruct *Crv2,
  51.                    CagdRType Blend)
  52. {
  53.     int i, j,
  54.     MaxAxis = CAGD_NUM_OF_PT_COORD(Crv1 -> PType),
  55.     Length = Crv1 -> Length,
  56.     Order = Crv1 -> Order;
  57.     CagdRType **NewPoints,
  58.     **Points1 = Crv1 -> Points,
  59.     **Points2 = Crv2 -> Points,
  60.     Blend1 = 1.0 - Blend;
  61.     CagdCrvStruct *NewCrv;
  62.  
  63.     if (Crv1 -> PType != Crv2 -> PType ||
  64.     Crv1 -> GType != Crv2 -> GType ||
  65.     Order != Crv2 -> Order ||
  66.     Length != Crv2 -> Length) {
  67.     SYMB_FATAL_ERROR(SYMB_ERR_CRVS_INCOMPATIBLE);
  68.     return NULL;
  69.     }
  70.     
  71.     NewCrv = CagdCrvNew(Crv1 -> GType, Crv1 -> PType, Length);
  72.     NewCrv -> Order = Order;
  73.     NewPoints = NewCrv -> Points;
  74.     if (Crv1 -> KnotVector != NULL)
  75.     NewCrv -> KnotVector = BspKnotCopy(Crv1 -> KnotVector, Length + Order);
  76.  
  77.     for (i = !CAGD_IS_RATIONAL_PT(Crv1 -> PType); i <= MaxAxis; i++) {
  78.     CagdRType
  79.         *Pts1 = &Points1[i][0],
  80.         *Pts2 = &Points2[i][0],
  81.         *NewPts = &NewPoints[i][0];
  82.  
  83.     for (j = Length - 1; j >= 0; j--)
  84.         *NewPts++ = *Pts1++ * Blend1 + *Pts2++ * Blend;
  85.     }
  86.  
  87.     return NewCrv;
  88. }
  89.  
  90. /*****************************************************************************
  91. * DESCRIPTION:                                                               M
  92. * Given two compatible curves (See function CagdMakeCrvsCompatible),         M
  93. * computes a morph between them using corner cutting approach.             M
  94. *   Returned is the new blended curve.                         M
  95. *                                                                            *
  96. * PARAMETERS:                                                                M
  97. *   Crv1, Crv2:  The two curves to blend.                                    M
  98. *   MinDist:     Minimal maximum distance between adjacent curves to make    M
  99. *                sure motion is visible. The curves will move at most twice  M
  100. *                that much in their maximal distance (roughly).             M
  101. *   SameLength:  If TRUE, length of curves is preserved, otherwise BBOX is.  M
  102. *   FilterTangenices: If TRUE, attempt is made to eliminate the intermediate M
  103. *                line representation.                         M
  104. *                                                                            *
  105. * RETURN VALUE:                                                              M
  106. *   CagdCrvStruct *:  The blended curve.                                     M
  107. *                                                                            *
  108. * KEYWORDS:                                                                  M
  109. *   SymbTwoCrvsMorphingCornerCut, morphing                            M
  110. *****************************************************************************/
  111. CagdCrvStruct *SymbTwoCrvsMorphingCornerCut(CagdCrvStruct *Crv1,
  112.                         CagdCrvStruct *Crv2,
  113.                         CagdRType MinDist,
  114.                         CagdBType SameLength,
  115.                         CagdBType FilterTangencies)
  116. {
  117.     int i, MorphListLength1, MorphListLength2,
  118.     Length = Crv1 -> Length;
  119.     CagdPointType
  120.         PType = Crv1 -> PType;
  121.     CagdPType Pl1, Pl2, Centroid1, Centroid2, TranslateFactor;
  122.     CagdVType Vl1, Vl2, VTmp;
  123.     CagdRType Scale1, Scale2, RotateFactor, *Centroid;
  124.     CagdCrvStruct *Crv, *Crv1MorphList, *Crv2MorphList;
  125.  
  126.     if (PType != Crv2 -> PType ||
  127.     Crv1 -> GType != Crv2 -> GType ||
  128.     Crv1 -> Order != Crv2 -> Order ||
  129.     Length != Crv2 -> Length) {
  130.     SYMB_FATAL_ERROR(SYMB_ERR_CRVS_INCOMPATIBLE);
  131.     return NULL;
  132.     }
  133.  
  134.     /* Make sure the curve is not closed. If so move first point epsilon. */
  135.     CagdCoerceToE3(Pl1, Crv1 -> Points, 0, PType);
  136.     CagdCoerceToE3(Pl2, Crv1 -> Points, Crv1 -> Length - 1, PType);
  137.     if (CGDistPointPoint(Pl1, Pl2) < LINE_EPSILON)
  138.         Crv1 -> Points[2][0] -= LINE_EPSILON * 10;
  139.     CagdCoerceToE3(Pl1, Crv2 -> Points, 0, PType);
  140.     CagdCoerceToE3(Pl2, Crv2 -> Points, Crv2 -> Length - 1, PType);
  141.     if (CGDistPointPoint(Pl1, Pl2) < LINE_EPSILON)
  142.         Crv2 -> Points[2][0] -= LINE_EPSILON * 10;
  143.  
  144.     /* Compute the morphing sequence to a line. */
  145.     Crv1MorphList = CrvMorphSeqCornerCut(Crv1, MinDist);
  146.     Crv2MorphList = CrvMorphSeqCornerCut(Crv2, MinDist);
  147.  
  148.     CagdCoerceToE3(Pl1, Crv1MorphList -> Points, 0, PType);
  149.     CagdCoerceToE3(Vl1, Crv1MorphList -> Points, Length - 1, PType);
  150.     PT_SUB(Vl1, Vl1, Pl1);
  151.     Scale1 = PT_LENGTH(Vl1);
  152.     PT_NORMALIZE(Vl1);
  153.  
  154.     CagdCoerceToE3(Pl2, Crv2MorphList -> Points, 0, PType);
  155.     CagdCoerceToE3(Vl2, Crv2MorphList -> Points, Length - 1, PType);
  156.     PT_SUB(Vl2, Vl2, Pl2);
  157.     Scale2 = PT_LENGTH(Vl2);
  158.     PT_NORMALIZE(Vl2);
  159.  
  160.     if (!SameLength) {
  161.     CagdRType Scale1a, Scale2a,
  162.         MidScale = sqrt(Scale1 * Scale2);
  163.     CagdBBoxStruct BBox1a, BBox1b, BBox2a, BBox2b;
  164.  
  165.     /* Preserve size as possible. */
  166.     CagdCrvBBox(Crv1MorphList, &BBox1a);
  167.     CagdCrvBBox((CagdCrvStruct *) CagdListLast(Crv1MorphList), &BBox1b);
  168.     CagdCrvBBox(Crv2MorphList, &BBox2a);
  169.     CagdCrvBBox((CagdCrvStruct *) CagdListLast(Crv2MorphList), &BBox2b);
  170.  
  171.     Scale1a = BBox1a.Max[0] - BBox1a.Min[0] >
  172.                   BBox1a.Max[1] - BBox1a.Min[1] ?
  173.         (BBox1b.Max[0] - BBox1b.Min[0]) / (BBox1a.Max[0] - BBox1a.Min[0]) :
  174.         (BBox1b.Max[1] - BBox1b.Min[1]) / (BBox1a.Max[1] - BBox1a.Min[1]);
  175.     Scale2a = BBox2a.Max[0] - BBox2a.Min[0] >
  176.                   BBox2a.Max[1] - BBox2a.Min[1] ?
  177.         (BBox2b.Max[0] - BBox2b.Min[0]) / (BBox2a.Max[0] - BBox2a.Min[0]) :
  178.         (BBox2b.Max[1] - BBox2b.Min[1]) / (BBox2a.Max[1] - BBox2a.Min[1]);
  179.  
  180.     MidScale *= sqrt(Scale1a * Scale2a);
  181.     Scale1 = MidScale / Scale1;
  182.     Scale2 = MidScale / Scale2;
  183.     }
  184.     else {
  185.     CagdRType
  186.         MidScale = sqrt(Scale1 * Scale2);
  187.  
  188.     /* Preserve length as possible. */
  189.  
  190.     Scale1 = MidScale / Scale1;
  191.     Scale2 = MidScale / Scale2;
  192.     }
  193.  
  194.     Centroid = ComputeCrvCentroid(Crv1);
  195.     PT_COPY(Centroid1, Centroid);
  196.     Centroid = ComputeCrvCentroid(Crv2);
  197.     PT_COPY(Centroid2, Centroid);
  198.     PT_SUB(TranslateFactor, Centroid2, Centroid1);
  199.  
  200.     CROSS_PROD(VTmp, Vl1, Vl2)
  201.     RotateFactor = atan2(VTmp[2], DOT_PROD(Vl1, Vl2));
  202.  
  203.     /* Apply the linear transform from one curve to the other, using the    */
  204.     /* Translation, scaling and rotation factors computed above.        */
  205.     MorphListLength1 = CagdListLength(Crv1MorphList);
  206.     for (i = 0, Crv = Crv1MorphList;
  207.      i < MorphListLength1;
  208.      i++, Crv = Crv -> Pnext) {
  209.     CagdRType
  210.         t = (MorphListLength1 - i) / ((CagdRType) MorphListLength1),
  211.         t2 = t / 2.0;
  212.     MatrixType RotateMat, ScaleMat, TranslateMat, Mat;
  213.  
  214.     Centroid = ComputeCrvCentroid(Crv);
  215.  
  216.     MatGenMatTrans(-Centroid[0], -Centroid[1], -Centroid[2], TranslateMat);
  217.     MatGenMatRotZ1(RotateFactor * t2, RotateMat);
  218.     MatMultTwo4by4(Mat, TranslateMat, RotateMat);
  219.  
  220.     MatGenMatUnifScale(1.0 + (Scale1 - 1.0) * t, ScaleMat);
  221.     MatMultTwo4by4(Mat, Mat, ScaleMat);
  222.  
  223.     MatGenMatTrans(Centroid1[0], Centroid1[1], Centroid1[2], TranslateMat);
  224.     MatMultTwo4by4(Mat, Mat, TranslateMat);
  225.  
  226.     MatGenMatTrans(TranslateFactor[0] * t2,
  227.                TranslateFactor[1] * t2,
  228.                TranslateFactor[2] * t2, TranslateMat);
  229.     MatMultTwo4by4(Mat, Mat, TranslateMat);
  230.      CagdCrvMatTransform(Crv, Mat);
  231.     }
  232.  
  233.     MorphListLength2 = CagdListLength(Crv2MorphList);
  234.     for (i = 0, Crv = Crv2MorphList;
  235.      i < MorphListLength2;
  236.      i++, Crv = Crv -> Pnext) {
  237.     CagdRType
  238.         t = (MorphListLength2 - i) / ((CagdRType) MorphListLength2),
  239.         t2 = t / 2.0;
  240.     MatrixType RotateMat, ScaleMat, TranslateMat, Mat;
  241.  
  242.     Centroid = ComputeCrvCentroid(Crv);
  243.  
  244.     MatGenMatTrans(-Centroid[0], -Centroid[1], -Centroid[2], TranslateMat);
  245.     MatGenMatRotZ1(-RotateFactor * t2, RotateMat);
  246.     MatMultTwo4by4(Mat, TranslateMat, RotateMat);
  247.  
  248.     MatGenMatUnifScale(1.0 + (Scale2 - 1.0) * t, ScaleMat);
  249.     MatMultTwo4by4(Mat, Mat, ScaleMat);
  250.  
  251.     MatGenMatTrans(Centroid2[0], Centroid2[1], Centroid2[2], TranslateMat);
  252.     MatMultTwo4by4(Mat, Mat, TranslateMat);
  253.  
  254.     MatGenMatTrans(-TranslateFactor[0] * t2,
  255.                -TranslateFactor[1] * t2,
  256.                -TranslateFactor[2] * t2, TranslateMat);
  257.     MatMultTwo4by4(Mat, Mat, TranslateMat);
  258.     CagdCrvMatTransform(Crv, Mat);
  259.     }
  260.  
  261.     if (FilterTangencies) {
  262.     TangencyFilterCrvMorph(&Crv1MorphList, &Crv2MorphList);
  263.     }
  264.  
  265.     /* Merge the two lists into one. */
  266.     Crv1MorphList = CagdListReverse(Crv1MorphList);
  267.     for (Crv = Crv1MorphList; Crv -> Pnext != NULL; Crv = Crv -> Pnext);
  268.     Crv -> Pnext = Crv2MorphList;
  269.  
  270.     DistanceFilterCrvMorph(Crv1MorphList, MinDist);
  271.  
  272.     return Crv1MorphList;
  273. }
  274.  
  275. /*****************************************************************************
  276. * DESCRIPTION:                                                               M
  277. * Given two compatible curves (See function CagdMakeCrvsCompatible),         M
  278. * computes a morph between them using multiresolution decomposition.         M
  279. *   Returned is a list of new blended curves.                     M
  280. *                                                                            *
  281. * PARAMETERS:                                                                M
  282. *   Crv1, Crv2:  The two curves to blend.                                    M
  283. *   BlendStep:   A step size of the blending.                                 M
  284. *                                                                            *
  285. * RETURN VALUE:                                                              M
  286. *   CagdCrvStruct *:  A list of blended curves.                              M
  287. *                                                                            *
  288. * KEYWORDS:                                                                  M
  289. *   SymbTwoCrvsMorphingMultiRes, morphing                               M
  290. *****************************************************************************/
  291. CagdCrvStruct *SymbTwoCrvsMorphingMultiRes(CagdCrvStruct *Crv1,
  292.                        CagdCrvStruct *Crv2,
  293.                        CagdRType BlendStep)
  294. {
  295.     CagdRType t;
  296.     CagdCrvStruct *MorphSeq, *TCrv;
  297.     SymbMultiResCrvStruct *CrvMultRes, *Crv1MultRes, *Crv2MultRes;
  298.  
  299.     if (Crv1 -> PType != Crv2 -> PType ||
  300.     Crv1 -> GType != Crv2 -> GType ||
  301.     Crv1 -> Order != Crv2 -> Order ||
  302.     Crv1 -> Length != Crv2 -> Length) {
  303.     SYMB_FATAL_ERROR(SYMB_ERR_CRVS_INCOMPATIBLE);
  304.     return NULL;
  305.     }
  306.     if (CAGD_IS_RATIONAL_CRV(Crv1)) {
  307.     SYMB_FATAL_ERROR(SYMB_ERR_RATIONAL_NO_SUPPORT);
  308.     return NULL;
  309.     }
  310.  
  311.     Crv1MultRes = SymbCrvMultiResDecomp(Crv1, FALSE),
  312.     Crv2MultRes = SymbCrvMultiResDecomp(Crv2, FALSE);
  313.     CrvMultRes = SymbCrvMultiResCopy(Crv1MultRes);
  314.     MorphSeq = CagdCrvCopy(Crv1);
  315.  
  316.     for (t = BlendStep; t < 1.0; t += BlendStep) {
  317.     int i;
  318.  
  319.     for (i = 0; i < Crv1MultRes -> Levels; i++) {
  320.         int l;
  321.         CagdCrvStruct
  322.         *Crv1MR = Crv1MultRes -> HieCrv[i],
  323.         *Crv2MR = Crv2MultRes -> HieCrv[i],
  324.         *CrvMR = CrvMultRes -> HieCrv[i];
  325.         CagdRType
  326.         **Points1 = Crv1MR -> Points,
  327.         **Points2 = Crv2MR -> Points,
  328.         **Points = CrvMR -> Points;
  329.  
  330.         for (l = 0; l < CrvMR -> Length; l++) {
  331.         InterpolateVector(Points1[1][l], Points1[2][l],
  332.                   Points2[1][l], Points2[2][l],
  333.                   &Points[1][l], &Points[2][l], t);
  334.         }
  335.     }
  336.  
  337.     TCrv = SymbCrvMultiResCompos(CrvMultRes);
  338.     LIST_PUSH(TCrv, MorphSeq);
  339.     }
  340.  
  341.     SymbCrvMultiResFree(CrvMultRes);
  342.     SymbCrvMultiResFree(Crv1MultRes);
  343.     SymbCrvMultiResFree(Crv2MultRes);
  344.  
  345.     return CagdListReverse(MorphSeq);
  346. }
  347. /*****************************************************************************
  348. * DESCRIPTION:                                                               *
  349. * Interpolate between vectors (X1, Y1) and (X2, Y2) by rotating and scaling. *
  350. *                                                                            *
  351. * PARAMETERS:                                                                *
  352. *   X1, Y1:   Coefficients of first vector.                                  *
  353. *   X2, Y2:   Coefficients of second vector.                                 *
  354. *   X, Y:     Place to substitute interpolated vector.                       *
  355. *   Blend:    Blending coefficients. Between zero and one.                   *
  356. *                                                                            *
  357. * RETURN VALUE:                                                              *
  358. *   void                                                                     *
  359. *****************************************************************************/
  360. static void InterpolateVector(CagdRType x1,
  361.                   CagdRType y1,
  362.                   CagdRType x2,
  363.                   CagdRType y2,
  364.                   CagdRType *x,
  365.                   CagdRType *y,
  366.                   CagdRType Blend)
  367. {
  368.     CagdRType Size1, Size2, Size, Angle, c, s;
  369.     CagdVType V1, V2, V;
  370.  
  371.     V1[0] = x1;
  372.     V1[1] = y1;
  373.     V1[2] = 0.0;
  374.     Size1 = PT_LENGTH(V1);
  375.     if (Size1 > PT_NORMALIZE_ZERO) {
  376.     V1[0] /= Size1;
  377.     V1[1] /= Size1;
  378.     }
  379.     V2[0] = x2;
  380.     V2[1] = y2;
  381.     V2[2] = 0.0;
  382.     Size2 = PT_LENGTH(V2);
  383.     if (Size2 > PT_NORMALIZE_ZERO) {
  384.     V2[0] /= Size2;
  385.     V2[1] /= Size2;
  386.     }
  387.  
  388.     CROSS_PROD(V, V1, V2)
  389.     Angle = atan2(V[2], DOT_PROD(V1, V2)) * Blend;
  390.     c = cos(Angle);
  391.     s = sin(Angle);
  392.     V[0] = V1[0] * c - V1[1] * s;
  393.     V[1] = V1[0] * s + V1[1] * c;
  394.     V[2] = 0.0;
  395.     Size = PT_LENGTH(V);
  396.     if (Size > PT_NORMALIZE_ZERO) {
  397.     V[0] /= Size;
  398.     V[1] /= Size;
  399.     }
  400.     *x = V[0] * (Size1 * (1.0 - Blend) + Size2 * Blend);
  401.     *y = V[1] * (Size1 * (1.0 - Blend) + Size2 * Blend);
  402. }
  403.  
  404. /*****************************************************************************
  405. * DESCRIPTION:                                                               *
  406. * Computes correspondance between the two sequences and eliminate all curves *
  407. * That a tangency correspondance is found.                     *
  408. *                                                                            *
  409. * PARAMETERS:                                                                *
  410. *   CrvSeq1, CrvSeq2:  Two curves morph3d into a line to filter out          *
  411. *                      tangencies.                                           *
  412. *                                                                            *
  413. * RETURN VALUE:                                                              *
  414. *   void                                                                     *
  415. *****************************************************************************/
  416. static void TangencyFilterCrvMorph(CagdCrvStruct **CrvSeq1,
  417.                    CagdCrvStruct **CrvSeq2)
  418. {
  419.     while ((*CrvSeq1) -> Pnext != NULL && (*CrvSeq2) -> Pnext != NULL) {
  420.     CagdRType *R;
  421.     CagdCrvStruct
  422.         *DCrv1 = CagdCrvDerive((*CrvSeq1) -> Pnext),
  423.         *DCrv2 = CagdCrvDerive((*CrvSeq2) -> Pnext),
  424.         *TanTan = SymbCrvDotProd(DCrv1, DCrv2);
  425.  
  426.     CagdCrvFree(DCrv1);
  427.     CagdCrvFree(DCrv2);
  428.  
  429.     R = SymbExtremumCntPtVals(TanTan -> Points, TanTan -> Length, TRUE);
  430.     if (R[1] > 0) {
  431.         CagdCrvStruct *TCrv;
  432.  
  433.         TCrv = *CrvSeq1;
  434.         *CrvSeq1 = (*CrvSeq1) -> Pnext;
  435.         CagdCrvFree(TCrv);
  436.         TCrv = *CrvSeq2;
  437.         *CrvSeq2 = (*CrvSeq2) -> Pnext;
  438.         CagdCrvFree(TCrv);
  439.     }
  440.     else
  441.         break;
  442.     }
  443. }
  444.  
  445. /*****************************************************************************
  446. * DESCRIPTION:                                                               *
  447. * Computes an estimate distance between two adjacent curves and equalize it. *
  448. *                                                                            *
  449. * PARAMETERS:                                                                *
  450. *   MorphSeq:  The sequence of morphed curves to filter to an approximate    *
  451. *              equal distance.                             *
  452. *   MinDist:   Distance between adjacent curves should be between this and   *
  453. *              twice this much.                             *
  454. *                                                                            *
  455. * RETURN VALUE:                                                              *
  456. *   void                                                                     *
  457. *****************************************************************************/
  458. static void DistanceFilterCrvMorph(CagdCrvStruct *MorphSeq,
  459.                    CagdRType MinDist)
  460. {
  461.     static CagdRType
  462.     Translate[3] = { 0.0, 0.0, 0.0 };
  463.     CagdCrvStruct *NextCrv;
  464.     CagdRType
  465.     IntMinDist = MinDist * MorphSeq -> Length;
  466.  
  467.     if (MorphSeq == NULL || MorphSeq -> Pnext == NULL)
  468.         return;
  469.  
  470.     while ((NextCrv = MorphSeq -> Pnext) != NULL) {
  471.     CagdRType IntDist;
  472.  
  473.     if ((IntDist = EstimateIntDistCrvCrv(MorphSeq, NextCrv)) <
  474.                                 IntMinDist &&
  475.         NextCrv -> Pnext != NULL) {
  476.         MorphSeq -> Pnext = NextCrv -> Pnext;
  477.         CagdCrvFree(NextCrv);
  478.     }
  479.     else if (IntDist > IntMinDist * 2.0) {
  480.         CagdCrvStruct
  481.             *MidCrv = SymbCrvAdd(MorphSeq, NextCrv);
  482.  
  483.         CagdCrvTransform(MidCrv, Translate, 0.5);
  484.  
  485.         MidCrv -> Pnext = NextCrv;
  486.         MorphSeq -> Pnext = MidCrv;
  487.     }
  488.     else
  489.         MorphSeq = MorphSeq -> Pnext;
  490.     }
  491. }
  492.  
  493.  
  494. /*****************************************************************************
  495. * DESCRIPTION:                                                               *
  496. * Computes an approximation to the given curve's centroid as an average of   *
  497. * its control points.                                                        *
  498. *                                                                            *
  499. * PARAMETERS:                                                                *
  500. *   Crv:        To compute centroid for.                                     *
  501. *                                                                            *
  502. * RETURN VALUE:                                                              *
  503. *   CagdRType *:   Centroid point.                                           *
  504. *****************************************************************************/
  505. static CagdRType *ComputeCrvCentroid(CagdCrvStruct *Crv)
  506. {
  507.     static CagdRType Centroid[3];
  508.     int i, j,
  509.     MaxAxis = CAGD_NUM_OF_PT_COORD(Crv -> PType),
  510.     Length = Crv -> Length;
  511.     CagdRType
  512.     **Points = Crv -> Points;
  513.  
  514.     for (j = 0; j < MaxAxis; j++)
  515.       Centroid[j] = 0.0;
  516.  
  517.     for (i = 0; i < Length; i++) {
  518.     for (j = 1; j <= MaxAxis; j++)
  519.         Centroid[j - 1] += Points[j][i];
  520.     }
  521.  
  522.     for (j = 0; j < MaxAxis; j++)
  523.     Centroid[j] /= Length;
  524.  
  525.     return Centroid;
  526. }
  527.  
  528.  
  529. /*****************************************************************************
  530. * DESCRIPTION:                                                               *
  531. * Computes a sequence of morphed curves from the given curve and upto a line *
  532. * using corner cutting. Curve is assumed to be at the origin.                *
  533. *                                                                            *
  534. * PARAMETERS:                                                                *
  535. *   Crv:        To corner cut to a line.                                     *
  536. *   MinDist:    Minimal maximum distance between adjacent curves to make     *
  537. *               sure motion is visible. The curves will move at most twice   *
  538. *               that much in their maximal distance (roughly).             *
  539. *                                                                            *
  540. * RETURN VALUE:                                                              *
  541. *   CagdCrvStruct *: List of curves representing the corner cut morph.         *
  542. *****************************************************************************/
  543. static CagdCrvStruct *CrvMorphSeqCornerCut(CagdCrvStruct *Crv,
  544.                        CagdRType MinDist)
  545. {
  546.     int i,
  547.     Length = Crv -> Length;
  548.     CagdPType Pt, Pl;
  549.     CagdVType Vl;
  550.     CagdCrvStruct
  551.     *LastCrv = CagdCrvCopy(Crv),
  552.     *MorphList = CagdCrvCopy(Crv);
  553.     CagdRType
  554.     IntMinDist = MinDist * Crv -> Length * 10.0;
  555.  
  556.     while (TRUE) {
  557.     CagdCrvStruct
  558.         *TCrv = CrvMorphOneCornerCut(LastCrv, MinDist);
  559.  
  560.     CagdCoerceToE3(Pl, TCrv -> Points, 0, TCrv -> PType);
  561.     CagdCoerceToE3(Vl, TCrv -> Points, TCrv -> Length - 1, TCrv -> PType);
  562.     PT_SUB(Vl, Vl, Pl);
  563.  
  564.     for (i = 1; i < Length - 1; i++) {
  565.         CagdCoerceToE3(Pt, TCrv -> Points, i, TCrv -> PType);
  566.         if (CGDistPointLine(Pt, Pl, Vl) > LINE_EPSILON)
  567.             break;
  568.     }
  569.     if (i >= Length - 1) {
  570.         LIST_PUSH(TCrv, MorphList);
  571.         break;
  572.     }
  573.     else {
  574.         CagdRType
  575.             DistBound = EstimateIntDistCrvCrv(TCrv, MorphList);
  576.  
  577.         CagdCrvFree(LastCrv);
  578.  
  579.         if (DistBound < IntMinDist) {
  580.         LastCrv = TCrv;
  581.         }
  582.         else {
  583.         LIST_PUSH(TCrv, MorphList);
  584.         LastCrv = CagdCrvCopy(TCrv);
  585.         }
  586.     }
  587.     }
  588.  
  589.     return MorphList;
  590. }
  591.  
  592. /*****************************************************************************
  593. * DESCRIPTION:                                                               *
  594. * Estimates the integrated distance between two curves as the sum of         *
  595. * distances between corresponding control points. Curves are assumed         *
  596. * compatible.                                     *
  597. *                                                                            *
  598. * PARAMETERS:                                                                *
  599. *   Crv1, Crv2:  Two curves to estimate integral distance between.           *
  600. *                                                                            *
  601. * RETURN VALUE:                                                              *
  602. *   CagdRType:   Integral distance bound.                                    *
  603. *****************************************************************************/
  604. static CagdRType EstimateIntDistCrvCrv(CagdCrvStruct *Crv1,
  605.                        CagdCrvStruct *Crv2)
  606. {
  607.     int i,
  608.     Length = Crv1 -> Length;
  609.     CagdRType
  610.     IntDist = 0.0,
  611.     **Points1 = Crv1 -> Points,
  612.     **Points2 = Crv2 -> Points;
  613.  
  614.     if (Crv1 -> PType != CAGD_PT_E2_TYPE && Crv1 -> PType != CAGD_PT_E3_TYPE) {
  615.     SYMB_FATAL_ERROR(SYMB_ERR_UNSUPPORT_PT);
  616.     return 0.0;
  617.     }
  618.     
  619.     for (i = 0; i < Length; i++)
  620.     {
  621.     IntDist += ABS(Points1[1][i] - Points2[1][i]) +
  622.                ABS(Points1[2][i] - Points2[2][i]) +
  623.            (Points1[3] != NULL ? ABS(Points1[3][i] - Points2[3][i])
  624.                        : 0.0);
  625.     }
  626.  
  627.     return IntDist;
  628. }
  629.  
  630. /*****************************************************************************
  631. * DESCRIPTION:                                                               *
  632. * Given a curve computes a smoother version of it using corner cutting.      *
  633. * The curve is assumed to be positioned around the origin.             *
  634. *                                                                            *
  635. * PARAMETERS:                                                                *
  636. *   Crv:     The curve to corner cut.                                    *
  637. *   Blend:       A parameter between zero and one. Zero will have zero       *
  638. *         affect on the curve, one the maximum corner cutting.         *
  639. *                                                                            *
  640. * RETURN VALUE:                                                              *
  641. *   CagdCrvStruct *:  The smoothed curve.                                    *
  642. *****************************************************************************/
  643. static CagdCrvStruct *CrvMorphOneCornerCut(CagdCrvStruct *Crv,
  644.                        CagdRType Blend)
  645. {
  646.     int i, j,
  647.     Length = Crv -> Length,
  648.     MaxAxis = CAGD_NUM_OF_PT_COORD(Crv -> PType);
  649.     CagdCrvStruct *SmoothedCrv;
  650.     CagdRType **Pts, **SmthPts, MinDist, TotalLength, TotalNewLength,
  651.         Translate[3],
  652.     *Nodes = CagdCrvNodes(Crv);
  653.  
  654.     SmoothedCrv = CagdCrvCopy(Crv);
  655.     SmthPts = SmoothedCrv -> Points;
  656.     Pts = Crv -> Points;
  657.  
  658.     MinDist = INFINITY;
  659.     TotalLength = 0.0;
  660.     for (i = 1; i < Length; i++) {
  661.     CagdRType
  662.         Dist = sqrt(SQR(SmthPts[1][i] - SmthPts[1][i - 1]) +
  663.             SQR(SmthPts[2][i] - SmthPts[2][i - 1]) +
  664.             SQR(SmthPts[3][i] - SmthPts[3][i - 1]));
  665.  
  666.     TotalLength += Dist;
  667.     if (MinDist > Dist)
  668.         MinDist = Dist;
  669.     }
  670.  
  671.     for (i = 1; i < Length - 1; i++) {
  672.     CagdRType
  673.         Dist1 = sqrt(SQR(Pts[1][i] - Pts[1][i - 1]) +
  674.              SQR(Pts[2][i] - Pts[2][i - 1]) +
  675.              SQR(Pts[3][i] - Pts[3][i - 1])),
  676.         Dist2 = sqrt(SQR(Pts[1][i] - Pts[1][i + 1]) +
  677.              SQR(Pts[2][i] - Pts[2][i + 1]) +
  678.              SQR(Pts[3][i] - Pts[3][i + 1])),
  679.         BlendNeighbor1 = Blend * SQR(MinDist / Dist2),
  680.         BlendNeighbor2 = Blend * SQR(MinDist / Dist1),
  681.         BlendSelf = 1.0 - BlendNeighbor1 - BlendNeighbor2;
  682.  
  683.     for (j = 1; j <= MaxAxis; j++)
  684.         SmthPts[j][i] = Pts[j][i] * BlendSelf +
  685.                     Pts[j][i - 1] * BlendNeighbor1 +
  686.                     Pts[j][i + 1] * BlendNeighbor2;
  687.     }
  688.     IritFree(Nodes);
  689.  
  690.     TotalNewLength = 0.0;
  691.     Translate[0] = Translate[1] = Translate[2] = 0.0;
  692.     for (i = 1; i < Length; i++) {
  693.     CagdRType
  694.         Dist = sqrt(SQR(SmthPts[1][i] - SmthPts[1][i - 1]) +
  695.             SQR(SmthPts[2][i] - SmthPts[2][i - 1]) +
  696.             SQR(SmthPts[3][i] - SmthPts[3][i - 1]));
  697.  
  698.     TotalNewLength += Dist;
  699.     }
  700.  
  701.     CagdCrvTransform(SmoothedCrv, Translate, TotalLength / TotalNewLength);
  702.  
  703.     return SmoothedCrv;
  704. }
  705.  
  706. /*****************************************************************************
  707. * DESCRIPTION:                                                               M
  708. * Given two compatible surfaces (See function CagdmakeSrfsCompatible),       M
  709. * computes a convex blend between them according to Blend which must be      M
  710. * between zero and one.                                 M
  711. *   Returned is the new blended surface.                     M
  712. *                                                                            *
  713. * PARAMETERS:                                                                M
  714. *   Srf1, Srf2:  The two surfaces to blend.                                  M
  715. *   Blend:       A parameter between zero and one                            M
  716. *                                                                            *
  717. * RETURN VALUE:                                                              M
  718. *   CagdSrfStruct *:   Srf2 * Blend + Srf1 * (1 - Blend).                    M
  719. *                                                                            *
  720. * KEYWORDS:                                                                  M
  721. *   SymbTwoSrfsMorphing, morphing                                            M
  722. *****************************************************************************/
  723. CagdSrfStruct *SymbTwoSrfsMorphing(CagdSrfStruct *Srf1,
  724.                    CagdSrfStruct *Srf2,
  725.                    CagdRType Blend)
  726. {
  727.     int i, j,
  728.     MaxAxis = CAGD_NUM_OF_PT_COORD(Srf1 -> PType),
  729.     ULength = Srf1 -> ULength,
  730.     VLength = Srf1 -> VLength,
  731.     UOrder = Srf1 -> UOrder,
  732.     VOrder = Srf1 -> VOrder;
  733.     CagdRType **NewPoints,
  734.     **Points1 = Srf1 -> Points,
  735.     **Points2 = Srf2 -> Points,
  736.     Blend1 = 1.0 - Blend;
  737.     CagdSrfStruct *NewSrf;
  738.  
  739.     if (Srf1 -> PType != Srf2 -> PType ||
  740.     Srf1 -> GType != Srf2 -> GType ||
  741.     UOrder != Srf2 -> UOrder ||
  742.     VOrder != Srf2 -> VOrder ||
  743.     ULength != Srf2 -> ULength ||
  744.     VLength != Srf2 -> VLength) {
  745.     SYMB_FATAL_ERROR(SYMB_ERR_SRFS_INCOMPATIBLE);
  746.     return NULL;
  747.     }
  748.     
  749.     NewSrf = CagdSrfNew(Srf1 -> GType, Srf1 -> PType, ULength, VLength);
  750.     NewSrf -> UOrder = UOrder;
  751.     NewSrf -> VOrder = VOrder;
  752.     NewPoints = NewSrf -> Points;
  753.     if (Srf1 -> UKnotVector != NULL)
  754.     NewSrf -> UKnotVector = BspKnotCopy(Srf1 -> UKnotVector,
  755.                         ULength + UOrder);
  756.     if (Srf1 -> VKnotVector != NULL)
  757.     NewSrf -> VKnotVector = BspKnotCopy(Srf1 -> VKnotVector,
  758.                         VLength + VOrder);
  759.  
  760.     for (i = !CAGD_IS_RATIONAL_PT(Srf1 -> PType); i <= MaxAxis; i++) {
  761.     CagdRType
  762.         *Pts1 = &Points1[i][0],
  763.         *Pts2 = &Points2[i][0],
  764.         *NewPts = &NewPoints[i][0];
  765.  
  766.     for (j = ULength * VLength - 1; j >= 0; j--)
  767.         *NewPts++ = *Pts1++ * Blend1 + *Pts2++ * Blend;
  768.     }
  769.  
  770.     return NewSrf;
  771. }
  772.